#include<iostream>
#include<ctime>
#include<cmath>
#include<fstream>
#include<algorithm>
using namespace std;

typedef struct node {
    char name[20];
    double x;
    double y;
}city;

static const int N = 30;//城市数量
static city citys[N];//城市列表
static double dic[N][N];//各城市之间距离
static int seq[N];//局部最优解记录的城市序列
static double answer;//局部最优解
const int tempterature = 1000;//初始温度
const double u = 0.998;//成功降温因子
const double v = 0.999;//失败降温因子
const int k = 100;//对每个温度迭代次数

double dic_two_point(city a, city b) {//根据给定两个城市的经纬度计算两城市间实际距离
    //    return sqrt(pow(a.x - b.x, 2) + pow(a.y - b.y, 2));
    const double PI = 3.1415926535898;
    double radLng1 = a.x * PI / 180.0;
    double radLng2 = b.x * PI / 180.0;
    double as = radLng1 - radLng2;
    double bs = (a.y-b.y) * PI / 180.0;
    double s = 2 * asin(sqrt(pow(sin(as/2), 2) + cos(radLng1) * cos(radLng2) * pow(sin(bs/2), 2)))* 6378.137;
    s =s * 10000 / 10000;
    return s;
}

void set_dic() {//计算各城市之间距离
    for (int i = 0; i<N; ++i) {
        for (int j = 0; j<N; ++j) {
            dic[i][j] = dic_two_point(citys[i], citys[j]);
        }
    }
}

double sum_dic(int* conf) {//计算路径总长度
    double temp = 0;
    for (int i = 1; i<N; ++i) {
        temp += dic_two_point(citys[conf[i]], citys[conf[i - 1]]);
    }
    temp += dic_two_point(citys[conf[0]], citys[conf[N - 1]]);
    return temp;
}

bool metropolis(double f1, double f2, double t) {//比较函数，是否用新解更新旧解
    if (f2 < f1)
        return true;
    double p = exp(-(f2 - f1) / t);
    int bignum = 1e9;
    if (rand() % bignum<p*bignum)
        return true;
    return false;
}

void generate(int* s) {//随机产生一组新解
    bool v[N];
    memset(v, false, sizeof(v));
    for (int i = 0; i<N; ++i) {
        s[i] = rand() % N;
        while (v[s[i]]) {
            s[i] = rand() % N;
        }
        v[s[i]] = true;
    }
}
void generateSwap(int* s) {//随机交换序列中的一组城市顺序
    int ti = rand() % N;
    int tj = ti;
    while (ti == tj)
        tj = rand() % N;
    for (int i = 0; i<N; ++i)
        s[i] = seq[i];
    swap(s[ti], s[tj]);
}

void  SAA() {//模拟退火算法（simulated annealing algorithm，SAA）实现
    double t = tempterature;
    int seq_t[N];
    for (int i = 0; i<N; ++i) {//初始化当前序列
        seq[i] = seq_t[i] = i;
    }
    double new_energy = 1, old_energy = 0;
    while (t>1e-9&&fabs(new_energy - old_energy)>1e-9) {//温度作为控制变量
        int t_k = k;
        int seq_tt[N];
        while (t_k--&&fabs(new_energy - old_energy)>1e-9) {//迭代次数作为控制变量
            generateSwap(seq_tt);
            new_energy = sum_dic(seq_tt);//new
            old_energy = sum_dic(seq_t);//old
            if (metropolis(old_energy, new_energy, t))
                for (int i = 0; i < N; ++i)
                    seq_t[i] = seq_tt[i];
        }
        new_energy = sum_dic(seq_t);//new
        old_energy = answer;//old
        if (metropolis(old_energy, new_energy, t)) {
            for (int i = 0; i < N; ++i)
                seq[i] = seq_t[i];
            answer = sum_dic(seq);
            t *= u;//接受新状态降温因子0.98
        }
        else
            t *= v;//不接受新状态降温因子0.99
    }
    answer = sum_dic(seq);
}
void input() {//读取城市经纬度数据
    ifstream ifile("ChinaCitys.txt");
    if (!ifile) {
        cout << "open field\n";
        return;
    }
    for(int i = 0; i < N; i++){
        ifile >> citys[i].name >> citys[i].x >> citys[i].y;
    }
}

void output() {//输出局部贪心最优解结果
    cout << "In this time, the best length is " << answer << ", the best road is : \n";
    for (int i = 0; i < N; ++i) {
        cout << citys[seq[i]].name<< " -> ";
    }
    cout << citys[seq[0]].name << endl;
}
int main() {
    clock_t start,finish;
    start=clock();
    srand(time(nullptr));
    input();//使用文件读取数据初始化城市经纬度
    set_dic();
    double best = DBL_MAX;
    for(int i = 0; i < 15; i++) {//依次从每个城市作为起始
        SAA();
        output();
        if(best > answer){
            best = answer;
        }
    }
    cout << "The current best length of the road is " << best << "km.\n";
    finish=clock();
    cout << "The run time is " << (double)(finish-start) * 1000 / CLOCKS_PER_SEC << "ms.\n";
    return 0;
}
